Все необходимые библиотеки:
library(readr) #task1
library(plotly) #task2
## Загрузка требуемого пакета: ggplot2
##
## Присоединяю пакет: 'plotly'
## Следующий объект скрыт от 'package:ggplot2':
##
## last_plot
## Следующий объект скрыт от 'package:stats':
##
## filter
## Следующий объект скрыт от 'package:graphics':
##
## layout
library(dplyr) #task3
##
## Присоединяю пакет: 'dplyr'
## Следующие объекты скрыты от 'package:stats':
##
## filter, lag
## Следующие объекты скрыты от 'package:base':
##
## intersect, setdiff, setequal, union
library(ggpubr)
library(rstatix)
##
## Присоединяю пакет: 'rstatix'
## Следующий объект скрыт от 'package:stats':
##
## filter
library(ggplot2) #task4
library(corrplot)
## corrplot 0.92 loaded
library(reshape2)
library(factoextra) #task5
## Welcome! Want to learn more? See two factoextra-related books at https://goo.gl/ve3WBa
library(pheatmap) #task6
library(FactoMineR) #task7
library(ggbiplot) #task8
## Загрузка требуемого пакета: plyr
## ------------------------------------------------------------------------------
## You have loaded plyr after dplyr - this is likely to cause problems.
## If you need functions from both plyr and dplyr, please load plyr first, then dplyr:
## library(plyr); library(dplyr)
## ------------------------------------------------------------------------------
##
## Присоединяю пакет: 'plyr'
## Следующие объекты скрыты от 'package:rstatix':
##
## desc, mutate
## Следующий объект скрыт от 'package:ggpubr':
##
## mutate
## Следующие объекты скрыты от 'package:dplyr':
##
## arrange, count, desc, failwith, id, mutate, rename, summarise,
## summarize
## Следующие объекты скрыты от 'package:plotly':
##
## arrange, mutate, rename, summarise
## Загрузка требуемого пакета: scales
##
## Присоединяю пакет: 'scales'
## Следующий объект скрыт от 'package:readr':
##
## col_factor
## Загрузка требуемого пакета: grid
№1 Загрузите датасет life_expectancy_data.RDS. Это данные с основными показателями, через которые высчитывается ожидаемая продолжительности жизни по метрике World Development Indicator на уровне стран. В данных оставлены строки, относящиеся к положению женщин в 2019 г.
# Чтение датасета
life_expectancy_data <- readRDS("life_expectancy_data.RDS")
# Просмотр первых строк датасета
head(life_expectancy_data)
## Country Year Gender Life expectancy Unemployment
## 1: Afghanistan 2019 Female 66.388 14.065000
## 2: Albania 2019 Female 80.201 11.322000
## 3: Algeria 2019 Female 78.133 18.629999
## 4: Angola 2019 Female 64.039 7.835000
## 5: Antigua and Barbuda 2019 Female 78.104 8.256285
## 6: Argentina 2019 Female 79.995 10.698000
## Infant Mortality GDP GNI
## 1: 42.9 18799450743 19098305631
## 2: 7.7 15400242875 15198661275
## 3: 18.6 171767403748 167568133680
## 4: 44.5 89417190341 81900890781
## 5: 5.1 1687533333 1581477852
## 6: 7.0 451932356086 434037404778
## Clean fuels and cooking technologies Per Capita
## 1: 36.0 494.1793
## 2: 80.7 5395.6595
## 3: 99.3 3989.6683
## 4: 49.6 2809.6261
## 5: 100.0 17376.6497
## 6: 99.8 10056.6379
## Mortality caused by road traffic injury Tuberculosis Incidence
## 1: 15.9 189
## 2: 11.7 16
## 3: 20.9 61
## 4: 26.1 351
## 5: 0.0 0
## 6: 14.1 29
## DPT Immunization HepB3 Immunization Measles Immunization Hospital beds
## 1: 66 66 64 0.4322222
## 2: 99 99 95 3.0523077
## 3: 91 91 80 1.8000000
## 4: 57 53 51 0.8000000
## 5: 95 99 93 2.5806250
## 6: 86 86 94 4.6100000
## Basic sanitation services Tuberculosis treatment Urban population
## 1: 49.00617 91.00000 25.754
## 2: 99.18307 88.00000 61.229
## 3: 86.13850 86.00000 73.189
## 4: 51.39362 69.00000 66.177
## 5: 85.52290 72.26667 24.506
## 6: 91.40773 47.00000 91.991
## Rural population Non-communicable Mortality Sucide Rate continent
## 1: 74.246 36.2 3.6 Asia
## 2: 38.771 6.0 2.7 Europe
## 3: 26.811 12.8 1.8 Africa
## 4: 33.823 19.4 2.3 Africa
## 5: 75.494 17.6 0.8 Americas
## 6: 8.009 12.1 3.3 Americas
Посмотрим структуру датасета:
str(life_expectancy_data)
## Classes 'data.table' and 'data.frame': 195 obs. of 23 variables:
## $ Country : chr "Afghanistan" "Albania" "Algeria" "Angola" ...
## $ Year : int 2019 2019 2019 2019 2019 2019 2019 2019 2019 2019 ...
## $ Gender : chr "Female" "Female" "Female" "Female" ...
## $ Life expectancy : num 66.4 80.2 78.1 64 78.1 ...
## $ Unemployment : num 14.06 11.32 18.63 7.84 8.26 ...
## $ Infant Mortality : num 42.9 7.7 18.6 44.5 5.1 ...
## $ GDP : num 1.88e+10 1.54e+10 1.72e+11 8.94e+10 1.69e+09 ...
## $ GNI : num 1.91e+10 1.52e+10 1.68e+11 8.19e+10 1.58e+09 ...
## $ Clean fuels and cooking technologies : num 36 80.7 99.3 49.6 100 ...
## $ Per Capita : num 494 5396 3990 2810 17377 ...
## $ Mortality caused by road traffic injury: num 15.9 11.7 20.9 26.1 0 ...
## $ Tuberculosis Incidence : num 189 16 61 351 0 29 26 2.2 6.9 6 ...
## $ DPT Immunization : num 66 99 91 57 95 ...
## $ HepB3 Immunization : num 66 99 91 53 99 ...
## $ Measles Immunization : num 64 95 80 51 93 ...
## $ Hospital beds : num 0.432 3.052 1.8 0.8 2.581 ...
## $ Basic sanitation services : num 49 99.2 86.1 51.4 85.5 ...
## $ Tuberculosis treatment : num 91 88 86 69 72.3 ...
## $ Urban population : num 25.8 61.2 73.2 66.2 24.5 ...
## $ Rural population : num 74.2 38.8 26.8 33.8 75.5 ...
## $ Non-communicable Mortality : num 36.2 6 12.8 19.4 17.6 ...
## $ Sucide Rate : num 3.6 2.7 1.8 2.3 0.8 ...
## $ continent : Factor w/ 5 levels "Africa","Americas",..: 3 4 1 1 2 2 4 2 5 4 ...
## - attr(*, ".internal.selfref")=<externalptr>
## - attr(*, "sorted")= chr "Country"
№2 Сделайте интерактивный plotly график любых двух нумерических колонок. Раскрасте по колонке континента, на котором расположена страна.
# Создание интерактивного графика с логарифмической осью Y и информационными подсказками
plot_ly(data = life_expectancy_data,
x = ~`Life expectancy`,
y = ~GDP,
type = 'scatter',
mode = 'markers',
color = ~continent,
colors = c("Africa" = 'blue', "Americas" = 'red', "Asia" = 'green', "Europe" = 'purple', "Oceania" = 'orange'),
marker = list(size = 10, opacity = 0.5),
text = ~paste('Country:', Country, '<br>GDP:', GDP, '<br>Life Expectancy:', `Life expectancy`),
# Добавляем подсказки
hoverinfo = 'text') %>%
layout(title = "Life Expectancy vs GDP by Continent (Log Scale)",
xaxis = list(title = "Life Expectancy"),
yaxis = list(title = "GDP (Log Scale)",
type = "log",
tickvals = c(1e5, 1e6, 1e7, 1e8, 1e9, 1e10, 1e11, 1e12, 1e13),
ticktext = c('100K', '1M', '10M', '100M', '1B', '10B', '100B', '1T', '10T'),
tickmode = 'array'),
legend = list(x = 1.05, y = 1, bordercolor = 'black', borderwidth = 1)) # Улучшаем выдачу при наведении
№3 Проведите тест на сравнение распределений колонки
Life expectancy между группами стран Африки и Америки. Вид
статистического теста определите самостоятельно. Визуализируйте
результат через библиотеку rstatix.
# Фильтрация данных для Африки и Америки
africa_america_data <- life_expectancy_data %>%
filter(continent %in% c("Africa", "Americas"))
# Проверка на нормальность распределений
africa_normality <- africa_america_data %>%
filter(continent == "Africa") %>%
pull(`Life expectancy`) %>%
shapiro_test()
america_normality <- africa_america_data %>%
filter(continent == "Americas") %>%
pull(`Life expectancy`) %>%
shapiro_test()
# Вывод результатов теста на нормальность
#print(africa_normality)
#print(america_normality)
# Выбор статистического теста в зависимости от нормальности распределения
if (africa_normality$p.value > 0.05 && america_normality$p.value > 0.05) {
# Если оба распределения нормальные, используем t-test
test_result <- t_test(`Life expectancy` ~ continent, data = africa_america_data)
} else {
# Если распределения не нормальные, используем тест Манна-Уитни
test_result <- wilcox_test(`Life expectancy` ~ continent, data = africa_america_data)
}
# Визуализация результатов
task_3 <- ggboxplot(africa_america_data, x = "continent", y = "Life expectancy",
color = "continent", palette = "jco") +
stat_compare_means(aes(label = ..p.signif..)) +
theme(legend.position = "none")
# Вывод графика
print(task_3)
№4 Сделайте новый датафрейм, в котором оставите все численные колонки
кромеYear. Сделайте корреляционный анализ этих данных.
Постройте два любых типа графиков для визуализации корреляций.
# Создание нового датафрейма
numeric_data <- life_expectancy_data %>%
select(-Year, -Country, -continent, -Gender) %>% # Исключаем нечисловые колонки и год
select_if(is.numeric)
# Вычисление корреляционной матрицы
cor_matrix <- cor(numeric_data, use = "pairwise.complete.obs") # 'pairwise.complete.obs' обрабатывает отсутствующие значения
# Визуализация корреляционной матрицы с помощью corrplot
corrplot(cor_matrix,
method = "circle",
type = "upper", # Показываем только верхнюю часть матрицы
order = "hclust", # Кластеризация переменных
tl.col = "black", # Цвет текста меток
tl.srt = 45, # Поворот текста меток
tl.cex = 0.6, # Размер текста меток
cl.cex = 0.8, # Размер текста цветовой шкалы
col = colorRampPalette(c("#BB0000", "#EE9999", "#FFFFFF", "#9999EE", "#0000BB"))(200),
addCoef.col = "black", # Цвет коэффициентов корреляции
addCoefasPercent = TRUE, # Показываем корреляцию в процентах
number.cex = 0.7, # Уменьшаем размер текста коэффициентов
diag = FALSE) # Убираем диагональные элементы
# Визуализация корреляционной матрицы с помощью ggplot2
library(reshape2) # Загрузка пакета reshape2 для преобразования матрицы
# Преобразование корреляционной матрицы для ggplot
cor_matrix_melt <- melt(cor_matrix)
names(cor_matrix_melt) <- c("Var1", "Var2", "value")
task_4 <- ggplot(cor_matrix_melt, aes(Var2, Var1, fill = value)) +
geom_tile(color = "white") +
scale_fill_gradient2(low = "blue", high = "red", mid = "white", midpoint = 0, limit = c(-1,1), space = "Lab") +
theme_minimal() +
theme(
axis.text.x = element_text(angle = 45, vjust = 1, size = 8, hjust = 1), # Поворачиваем текст на оси X и уменьшаем размер
axis.text.y = element_text(size = 8), # Уменьшаем размер текста на оси Y
axis.title.x = element_blank(), # Убираем подпись оси X
axis.title.y = element_blank(), # Убираем подпись оси Y
legend.title = element_text(size = 10), # Уменьшаем размер заголовка легенды
legend.text = element_text(size = 9) # Уменьшаем размер текста легенды
) +
labs(
fill = "Корреляция\nПирсона" # Задаем заголовок легенды
) +
coord_fixed()
print(task_4)
№5 Постройте иерархическую кластеризацию на этом датафрейме.
# Подготовка данных для кластеризации
numeric_data <- life_expectancy_data %>%
select(-Year, -Country, -continent, -Gender) %>%
select_if(is.numeric)
# Создаём матрицу дистанций для числовых данных
numeric_data_dist <- dist(numeric_data, method = "euclidean")
# Выводим часть матрицы дистанций
as.matrix(numeric_data_dist)[1:6, 1:6]
## 1 2 3 4 5 6
## 1 0 5173184728 213172429004 94504126087 24487894415 599815446359
## 2 5173184728 0 218328067643 99637823809 19325270693 604967915615
## 3 213172429004 218328067643 0 118829433029 237652545186 386650065250
## 4 94504126087 99637823809 118829433029 0 118944107933 505388335722
## 5 24487894415 19325270693 237652545186 118944107933 0 624290420516
## 6 599815446359 604967915615 386650065250 505388335722 624290420516 0
# Вычисляем дендрограмму кластеров
numeric_data_hc <- hclust(d = numeric_data_dist, method = "ward.D2")
# Визуализируем дендрограмму
fviz_dend(numeric_data_hc, cex = 0.1) # Размер меток на дендрограмме
## Warning: The `<scale>` argument of `guides()` cannot be `FALSE`. Use "none" instead as
## of ggplot2 3.3.4.
## ℹ The deprecated feature was likely used in the factoextra package.
## Please report the issue at <https://github.com/kassambara/factoextra/issues>.
## This warning is displayed once every 8 hours.
## Call `lifecycle::last_lifecycle_warnings()` to see where this warning was
## generated.
№6 Сделайте одновременный график heatmap и иерархической кластеризации.
Содержательно интерпретируйте результат.
pheatmap(numeric_data,
show_rownames = FALSE,
clustering_distance_rows = "euclidean",
clustering_distance_cols = "euclidean",
clustering_method = "ward.D2",
cutree_rows = 5,
cutree_cols = length(colnames(numeric_data)),
angle_col = 45,
main = "Dendrograms for clustering rows and columns with heatmap")
Построена heatmap с дендрограммами для строк и столбцов на основе всех
количественных переменных датасета.
Основная проблема на мой взгляд - масштаб значений. Цветовая шкала на графике показывает очень большой диапазон значений, необходима нормализация или стандартизация данных перед их визуализацией. Без этого шага переменные с большими значениями, такие как ВВП, доминируют в анализе и искажают картину кластеризации.
Иерархическая кластеризация была выполнена как для строк, так и для столбцов. Дендрограмма слева показывает группировку переменных, которые похожи по своим паттернам в разных странах или регионах.
Переменные ВВП (GDP) и уровень жизни (GNI), как правило, имеют высокие значения, что может указывать на их сильное влияние на другие измерения здоровья и благополучия.
Попробем улучшить график:
# Нормализация данных
numeric_data_normalized <- as.data.frame(lapply(numeric_data, function(x) {
(x - min(x)) / (max(x) - min(x))
}))
# Визуализация стандартизированных нормализованных данных
pheatmap(numeric_data_normalized,
show_rownames = FALSE,
clustering_method = "ward.D2",
angle_col = 45,
main = "Dendrograms with Heatmap (Normalized Data)")
Не уверена, что нормализация - лучший выход для этих данных, но теперь
можно анализировать.
На графике видно, что экономические показатели, такие как ВВП (GDP - Gross Domestic Product), уровень жизни (GNI - Gross National Income) и доход на душу населения (Per Capita), образуют отдельный кластер. Это указывает на то, что эти переменные, как правило, коррелируют между собой.
№7 Проведите PCA анализ на этих данных. Проинтерпретируйте результат.
# Выполнение PCA на числовых данных с стандартизацией
life_expectancy_pca <- prcomp(numeric_data, scale = TRUE)
# Оценка результатов PCA
summary(life_expectancy_pca)
## Importance of components:
## PC1 PC2 PC3 PC4 PC5 PC6 PC7
## Standard deviation 2.7526 1.4841 1.3952 1.17177 1.08375 0.96347 0.9288
## Proportion of Variance 0.3988 0.1159 0.1025 0.07227 0.06182 0.04886 0.0454
## Cumulative Proportion 0.3988 0.5147 0.6172 0.68945 0.75126 0.80012 0.8455
## PC8 PC9 PC10 PC11 PC12 PC13 PC14
## Standard deviation 0.85740 0.69263 0.68937 0.59106 0.54986 0.47085 0.36596
## Proportion of Variance 0.03869 0.02525 0.02501 0.01839 0.01591 0.01167 0.00705
## Cumulative Proportion 0.88421 0.90946 0.93447 0.95286 0.96877 0.98044 0.98749
## PC15 PC16 PC17 PC18 PC19
## Standard deviation 0.34546 0.26941 0.20224 0.06968 1.017e-15
## Proportion of Variance 0.00628 0.00382 0.00215 0.00026 0.000e+00
## Cumulative Proportion 0.99377 0.99759 0.99974 1.00000 1.000e+00
# Визуализация доли объясненной дисперсии для каждой главной компоненты
fviz_eig(life_expectancy_pca, addlabels = TRUE, ylim = c(0, 40))
Результаты PCA показывают, что первые четыре главные компоненты (PC1,
PC2, PC3, PC4) объясняют примерно 68.945% общей вариации в данных.
Первая компонента (PC1) объясняет почти 40% вариации, последующие
компоненты (PC2, PC3, PC4) добавляют ещё около 29% к объяснению
вариации.
# Визуализация вклада переменных
fviz_pca_var(life_expectancy_pca, col.var = "contrib")
По этому графику можно выделить группы переменных, но лучше корректно их визуализируем в следующем задании)
№8 Постройте biplot график для PCA. Раскрасьте его по значениям
континентов.Переведите его в plotly. Желательно, чтобы при
наведении на точку, вы могли видеть название страны.
дальше будет несколько попыток совместить все тебования в одном графике, пока не получилось совместить
# Строим biplot с помощью ggbiplot
ggbiplot(life_expectancy_pca,
scale = 0, # Указываем scale = 0, чтобы использовать оригинальные масштабы
alpha = 0.1 # Установка прозрачности точек
) +
theme_minimal() # Применение минималистичной темы
Здесь не получилось изменить numeric формат континентов из-за разницы длинны векторов при таком подходе:
data <- life_expectancy_data %>%
filter(!is.na(continent) &
!is.na(`Life expectancy`) &
`Life expectancy` != 0 &
Unemployment != 0 &
`Infant Mortality` != 0 &
GDP != 0 &
GNI != 0 &
`Clean fuels and cooking technologies` != 0 &
`Per Capita` != 0 &
`Mortality caused by road traffic injury` != 0 &
`Tuberculosis Incidence` != 0 &
`DPT Immunization` != 0 &
`HepB3 Immunization` != 0 &
`Measles Immunization` != 0 &
`Hospital beds` != 0 &
`Basic sanitation services` != 0 &
`Tuberculosis treatment` != 0 &
`Urban population` != 0 &
`Rural population` != 0 &
`Non-communicable Mortality` != 0 &
`Sucide Rate` != 0
)
# Преобразование всех столбцов в числовой формат
data <- as.data.frame(sapply(data, as.numeric))
## Warning in lapply(X = X, FUN = FUN, ...): в результате преобразования созданы
## NA
## Warning in lapply(X = X, FUN = FUN, ...): в результате преобразования созданы
## NA
# PCA
pca <- prcomp(data[, -(1:3)], scale = TRUE)
ggbiplot(pca,
scale = 0,
groups = data$continent, # Раскраска по континентам
alpha = 0.2) +
theme_minimal()
Меняем подход:
# Выбераем только числовые переменные для анализа PCA
numeric_data <- life_expectancy_data[, c("Life expectancy", "Unemployment",
"Infant Mortality", "GDP","GNI",
"Clean fuels and cooking technologies", "Per Capita",
"Mortality caused by road traffic injury", "Tuberculosis Incidence",
"DPT Immunization", "HepB3 Immunization", "Measles Immunization",
"Hospital beds", "Basic sanitation services", "Tuberculosis treatment",
"Urban population", "Rural population", "Non-communicable Mortality",
"Sucide Rate")]
# PCA
pca_result <- prcomp(numeric_data, scale. = TRUE)
# Дадафрейм данных с результатами PCA и названиями стран
pca_df <- data.frame(PC1 = pca_result$x[, 1], PC2 = pca_result$x[, 2], Country = life_expectancy_data$Country)
# Добавляем информацию о континентах к датафрейму данных
pca_df <- cbind(pca_df, Continent = life_expectancy_data$continent)
# Конвертируем в plotly
biplotly <- plot_ly(pca_df, x = ~PC1, y = ~PC2, text = ~Country, color = ~Continent) %>%
add_markers() %>%
layout(title = "Biplot график для PCA с раскраской по континентам",
xaxis = list(title = "Главная компонента 1"),
yaxis = list(title = "Главная компонента 2"),
hovermode = "closest")
biplotly
№10. Сравните результаты отображения точек между алгоритмами PCA и UMAP.
Так бы я вывела 2 графика для сравнения, если бы смогла решить проблему с continent:
library(umap)
#UMAP
umap_result <- umap(numeric_data, n_neighbors = 10, n_components = 2, metric = "euclidean")
# Датафрейм с результатами UMAP и названиями стран
umap_df <- data.frame(UMAP1 = umap_result$layout[, 1], UMAP2 = umap_result$layout[, 2], Country = life_expectancy_data$Country)
# biplot для PCA
pca_plot <- ggplot(pca_df, aes(x = PC1, y = PC2, text = Country, color = Continent)) +
geom_point() +
labs(title = "Biplot график для PCA",
x = "Главная компонента 1", y = "Главная компонента 2") +
theme_minimal()
# scatter plot для UMAP
umap_plot <- ggplot(umap_df, aes(x = UMAP1, y = UMAP2, text = Country, color = Continent)) +
geom_point() +
labs(title = "Scatter plot для UMAP",
x = "UMAP 1", y = "UMAP 2") +
theme_minimal()
# ggplotly
pca_plotly <- ggplotly(pca_plot)
umap_plotly <- ggplotly(umap_plot)
pca_plotly
umap_plotly